program analysis
implicit none
integer::i,j,k
integer::betanum
integer,parameter::nbeta=48,nlatt=1728
integer,parameter::nrec=1
integer,parameter::nswp1=2000,nexch1=200
integer,parameter::nswp2=0,nexch2=500
double precision,dimension(0:nbeta-1)::beta
integer(kind=4),parameter::nsample=nexch1*nswp1+nexch2*nswp2
integer(kind=4),parameter::nsample1=nsample
integer(kind=4),parameter::nsample2=nsample
integer,parameter::selectflag=0
integer(kind=4),dimension(0:nbeta-1)::samcount
double precision,dimension(0:nbeta-1)::e_ave,u_ave,ee_ave,uu_ave,spheat_ave,psus_ave,m_ave
double precision,dimension(0:nbeta-1)::absu_ave,psusm_ave,mm_ave,chi
double precision,dimension(0:nbeta-1)::proe,prou,prom
integer,parameter::uflag=0

 character(len=3)::filename
 integer::filenum1,filenum2,filenum3
 integer::my_rank,ifile


 do my_rank=0,nbeta-1
 filenum1=my_rank/100
 filenum2=(my_rank-filenum1*100)/10
 filenum3=(my_rank-filenum1*100-filenum2*10)
 filename=achar(filenum1+48)//achar(filenum2+48)//achar(filenum3+48)
 open(my_rank+1000,file='../outdata/'//filename//'.bin',form='unformatted')
 end do
  write(*,*)'the files are opened corretly'

samcount(:)=1
e_ave(:)=0.0
u_ave(:)=0.0
m_ave(:)=0.0
mm_ave(:)=0.0
chi(:)=0.0
ee_ave(:)=0.0
uu_ave(:)=0.0
absu_ave(:)=0.0
psusm_ave(:)=0.0

open(200,file='results.dat')
!read the initial beta file
  open(600,file='../beta_init.dat')
  do i=0,nbeta-1
  read(600,*)k,beta(i)
  end do
  close(600)
  write(*,*)'read the initial beta is ok'

!read the data of replica_exchange from the binary file
do ifile=1000,1000+nbeta-1
    
    if(selectflag==1)then
    do j=1,nsample1
       read(ifile)betanum
    end do
    end if
    
    do i=1,nsample2
       if(uflag==1)then
       read(ifile)betanum,proe(betanum),prou(betanum)
       else
       read(ifile)betanum,proe(betanum),prom(betanum)
       end if
       e_ave(betanum)=e_ave(betanum)+proe(betanum)
       u_ave(betanum)=u_ave(betanum)+prou(betanum)
       m_ave(betanum)=m_ave(betanum) + prom(betanum)
       absu_ave(betanum)=absu_ave(betanum)+abs(prou(betanum))
       ee_ave(betanum)=ee_ave(betanum)+proe(betanum)**2.0
       uu_ave(betanum)=uu_ave(betanum)+prou(betanum)**2.0        
       mm_ave(betanum)=mm_ave(betanum)+prom(betanum)**2.0
       if(samcount(betanum)>nsample2)then
       write(*,*)'the number is error'
       call abort
       end if
       samcount(betanum)=samcount(betanum)+1
    end do
    write(*,*)'reading the',ifile,'is ok' 
end do

do i=0,nbeta-1
 if(samcount(i)/=nsample2+1)then
 write(*,*)'samcount is error'
 call abort
 end if
end do
write(*,*)'reading the data is ok'

open(800,file='./energy.dat')
open(801,file='./spheat.dat')
open(802,file='./magnetism.dat')
open(803,file='./susceptibility.dat')
do i=0,nbeta-1
  e_ave(i)=e_ave(i)/nsample2
  u_ave(i)=u_ave(i)/nsample2
  m_ave(i)=m_ave(i)/nsample2
  mm_ave(i)=mm_ave(i)/nsample2
  ee_ave(i)=ee_ave(i)/nsample2
  uu_ave(i)=uu_ave(i)/nsample2
  absu_ave(i)=absu_ave(i)/nsample2
  spheat_ave(i)=(ee_ave(i)-e_ave(i)**2.0)*(beta(i)**2.0)
  chi(i)=(mm_ave(i)- m_ave(i)**2.0)*(beta(i)**1.0)
  psus_ave(i)=(uu_ave(i)-u_ave(i)**2.0)*beta(i)
  psusm_ave(i)=(uu_ave(i)-absu_ave(i)**2.0)*beta(i)
  write(200,'(5f20.8)') 1/beta(i),e_ave(i)/nlatt,spheat_ave(i)/nlatt,m_ave(i)/nlatt,chi(i)/nlatt
  write(800, '(2F20.8)') 11.59/beta(i), e_ave(i)/nlatt
  write(801, '(2F20.8)') 11.59/beta(i), spheat_ave(i)/nlatt
  write(802, '(2F20.8)') 11.59/beta(i), m_ave(i)/nlatt
  write(803, '(2F20.8)') 11.59/beta(i), chi(i)/nlatt

end do

close(200)
close(800)
close(801)
close(802)
close(803)
end program
